// Problem:    153
// Title:      Investigating Gaussian Integers.
// Date:       05-May-2007
// Difficulty: 2
// Fun Rating: 1
// Time:       n
// Space:      n^0.5
// Closed:     Y
// Answer:     17971254122360635

#define ENABLE_CACHE    1
#define ENABLE_CALLSTAT 0

#include <iostream>
#include <cmath>
#if ENABLE_CALLSTAT
#include <map>
#endif
#if ENABLE_CACHE
#include <vector>
#endif

#if ENABLE_CALLSTAT
static std::map<int,int> call_stat;
#endif

#if ENABLE_CACHE
static std::vector<long long> R_cache;
#endif

// Computes the sum of real divisors of all integers n <= N.
static long long R(int N)
{
#if ENABLE_CALLSTAT
	++call_stat[N];
#endif

#if ENABLE_CACHE
	// Try find the value from cache.
	if (N < (int)R_cache.size() && R_cache[N] > 0)
		return R_cache[N];
#endif

	// Computes R(N).
	long long sum = 0;
	for (int a = (int)sqrt((double)N); a >= 1; a--)
	{
		int m = N / a;
		long long s = (long long)(3*a+m+1)*(m-a)/2 + a;
		sum += s;
	}

#if ENABLE_CACHE
	// Update cache.
	if (N < (int)R_cache.size())
		R_cache[N] = sum;
#endif

	return sum;
}

// Computes the sum of complex divisors of integer n <= N.
static long long C(int N)
{
	// Initilize the sum as generated by (1,1).
	long long sum = 2*R(N/2);

	// Initialize a farey sequence of order sqrt(N-1).
	int r = (int)sqrt((double)(N-1));
	int pp = 0, qq = 1;
	int p = 1, q = r;

	while (q > 1)
	{
		// Compute the next term in the farey sequence.
		int t = (r+qq)/q;
		int p_next = t*p-pp;
		int q_next = t*q-qq;
		
		// Compute the sum of divisors generated from this (p, q).
		int m = N/(p*p+q*q);
		if (m > 0)
		{
			sum += 2*(p+q)*R(m);
		}
		else
		{
			r--;
		}

		pp = p;
		qq = q;
		p = p_next;
		q = q_next;
	}
	return sum;
}

#define USE_NEW 1

void solve_problem_153()
{
#if 1
	const int N = 100000000;
#else
	const int N = 100000;
	//const int N = 5;
#endif

	// Initialize R_cache
	R_cache.resize((int)sqrt((double)N)+1);

	// Compute the sum of real divisors.
	long long s_real = R(N);
	//std::cout << "R(N) = " << s_real << std::endl;

	// Compute the sum of complex divisors.
	long long s_complex = C(N);
	//std::cout << "C(N) = " << s_complex << std::endl;

	// Display result.
	std::cout << (s_real+s_complex) << std::endl;

#if ENABLE_CALLSTAT
	// Display call stat
	int max_display = 200;
	for (auto it = call_stat.cbegin(); it != call_stat.cend(); ++it)
	{
		if (--max_display < 0)
			break;
		std::cout << it->first << " " << it->second << std::endl;
	}
#endif
}
